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Abstract 

An inhomogeneous random recursive lattice is constructed from the multi- 
branched Husimi square lattice. The number of repeating units connected 
on one vertex is randomly set to be 2 or 3 with a fixed ratio P 2 or P3 with 
P 2 + P3 = 1. The lattice is designed to describe complex thermodynamic 
systems with variable coordinating neighbors, e.g. the asymmetric range 
around the surface of a bulk system. Classical ferromagnetic spin-1 Ising 
model is solved on the lattice to achieve an annealed solution via the local 
exact calculation technique. The model exhibits distinct spontaneous mag¬ 
netization similar to the deterministic system, with however rigorous thermal 
fluctuations and significant singularities on the entropy behavior around the 
critical temperature, indicating a complex superheating frustration in the 
cross-dimensional range induced by the stochasticity. The critical tempera¬ 
ture was found to be exponentially correlated to the structural ratio P with 
the coefficient fitted as 0.53187, while the ground state energy presents linear 
correlation to P, implying a well-defined average property according to the 
structural ratio. 
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1. Introduction 


The Bethe or Bethe-like recursive lattices generally refer to the fractal ar¬ 
rangements of repeating units recursively connected to neighbors only on the 
sharing vertex, with no connection bond lies crossing layers. It has become a 
powerful methodology in various fields such as thermodynamics [2, 2], graph 
theory^, optimization problems | lj. [•'>[ and so on. In statistical physics, one 
important application of the recursive lattice is to approximate the regular 
lattice with the identical coordination number to solve a thermodynamic sys¬ 
tem (e.g. Ising model) on it. As one of few exactly calculable models, it has 
been proven to be a reliable method |6|, with the advantage of exact calcu¬ 
lation and simple iterative approach [7j, to be applied in numerous physical 
systems, e.g. alloy j^|, spin glass jo[, polymers jlo[, biomacromolecule H3 
etc. 


As a versatile extension of the Bethe lattice assembled by single dots 
and bonds, the Husirni lattice employing simple geograp hic shapes, such as 


square, triangle, tetrahedron, hexagon, or cube [12Hl4j has also been de 


veloped for decades to describe various systems with particular geographic 
properties |15|-il7|. Similar to the Bethe lattice, The independence feature of 


units enables the exact calculation on Husirni lattice regardless of the dimen¬ 
sions of the geographic unit, and mean field approximation is unnecessary 
since the interactions are confined within a unit and not shared by others. 
The calculation method usually relies on the recursive approach, which is 
featured as simple and less computation effort costing thanks to the homo¬ 


geneous self-similar structure [13|, [18 


Nevertheless, the recursive feature is accompanied with several disadvan¬ 
tages of this lattice methodology. Firstly, the repeating structure implies a 
homogeneous system, it is only suitable to describe systems of uniform tex¬ 
ture. Some particular however important cases, e.g. the confined geometry 
or structural transformation, are enormously difficult, if not impossible, to 
be simulated by the recursive lattices. Therefore, besides a few investiga¬ 
tions on the thermodynamics on the surface/thin film employed moderately 


inhomogeneous structure to present the boundary of a bulk system [19|, [20 


the reports on the application of recursive lattice onto inhomogeneous sys¬ 
tems were very rare. Secondly, recursive lattice is considered to be a reliable 
approximation to regular lattice based on the identical coordination num¬ 
ber q. Therefore, the manipulation of coordination number(s) is critical in 
constructing a recursive lattice for particular requirements. While it is easy 
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to draw a regular lattice with an arbitrary q, achieving an odd q in recur¬ 
sive lattice usually requires awkward design of unit selection and branch 
number, and even worse a prime number of q is impossible. Furthermore, 
when the randomness is necessary in a recursive lattice model, the common 
method is to add random terms in the Hamiltonian, e.g. a random external 
field as noise or random exchange couplings parameter j 2 l|, 22 ], while the 


structural randomness is difficult to be presented due to the homogeneity of 
recursions. 

Therefore, it is expected to be a considerable contribution to this held 
that if new designs of recursive lattice and calculation methods are developed 
to handle the above concerns, and make the recursive methodology more ver¬ 
satile in describing inhomogeneous systems. Recently we reported an Husimi 
lattice of random square-cube recursion, on which the simple Ising model can 
be solved by conventional exact calculation technique with moderate modifi¬ 
cation and exhibits well-defined thermal behaviors j23|. Following the same 
principle, in this work we have developed a randomly multi-branched Husimi 
square lattice and solved the simple ferromagnetic Ising model on it. The 
lattice is featured by randomly two or three square units connected on one 
vertex, then an inhomogeneous system of variable coordination numbers can 
be achieved with the identical unit cells. The spontaneous magnetization 
with critical temperature Tq, thermodynamics around the singularity, and 
the thermal fluctuation caused by stochastic structure, were investigated with 
the variation of structural ratio. 


2. Modeling and Calculation 

2.1. Lattice Construction 

The original Husimi lattice was a tree-like graph assembled by squares 
with two units connecting on one vertex. Since its development, derivative 
structures of three and more squares connection have been also investigated. 
In this way, to achieve an inhomogeneous structure with variable q, it is a 
natural choice to have random number of branches connected in the lattice. 
To keep the investigation simple, in this work we will only study the randomly 
2 or 3-branched Husimi square lattice with q = 4 or 6 as demonstrated in 
Fig{T| A structural ratio P 2 or P 3 can be defined to indicate the probability 
to have 2 or 3 branches on one vertex, with obviously P2 + P3 — 1. However, 
as a probable reason why this type of lattice has not been reported before, 
the random structure destroys the recursive homogeneity and then makes the 
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Figure 1: A demonstration of Husirni lattice with randomly two or three 
branches joint on one vertex. 

iterative approach unfeasible, therefore the lattice shown in FigJT] is not the 
actual model studied in this paper, and particular limitations on the structure 
must be applied to achieve an exact calculation, which will be detailed later. 

While the branch number is randomly 2 or 3 in this lattice, with the 
structural ratio P we can define an “analog branch number” as 


( 1 ) 


analog B = 2xPo + 3x P 3 . 


to present the average branch number of the lattice, and similarly an analog 
coordination number q can also be defined as analog q = 2 x analog B. By 
this means, the term “cross-dimensional” in this paper refers to the gradient 
of analog q: Taking a finite regular cubic lattice with a 2D surface as an 
example, we have q — 5 on the surface and q — 6 in the bulk, then in the 
near-surface region a randomly sampled site will has a probability to be of 
either q depending on the depth. Therefore, a gradually variation of analog 
q from 5 to 6 well represents the cross-dimensional range from the surface 
to bulk in this case. Similarly, the case of a 2D layer crossing to thin him 
can be described by the variation of analog q from 4 to 6 . For an additional 
clarification, both 3-branched Husirni square lattice and Husirni cubic lattice 
have been proved to be a good approximation to the regular cubic lattice 
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The simplest ferromagnetic spinil Ising model was applied on the lattice 
in this paper: 



( 2 ) 


<i,j> 


4 




without external magnetic field H. The weights of one configuration 7 of a 
square unit is given by 
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w( 7 ) = exp(-/3 ^ - JSiSj ), 


( 3 ) 


<i,j> 


where /3 is the reversal temperature as l/Zc^T, the Boltzmann constant ke 
is set as one. We have the partition function of the entire system as 


z = 


( 4 ) 


r a 


where the T = (££) 7 a denotes the state of the lattice as an ensemble of unit 
a. 

In this paper we setup a uniform ferromagnetic coupling J %] = 1, then the 
state of system only depends on the structural properties. Without external 
magnetic field, we can expect a half-half probability of spin state on each site 
at high temperature, a uniform orientation pointing to either up or down of 
all spins at low temperature, and a spontaneous magnetization occurring in 
between. The only question being focused on in this paper is that, how this 
transition behaves in the cross-dimensional situation on an inhomogeneous 
lattice. 

2.2. Partial Partition Function and Cavity Field 

The lattice is designed of infinite size, nevertheless for an iterative ap¬ 
proach it is necessary to imagine an original point where the entire lattice 
contribute to. Furthermore, the structure must be symmetrical to the origi¬ 
nal point, and subsequently the symmetry of sub-trees contributing onto one 
unit is required, otherwise the unique structure of an arbitrary sub-tree is 
impossible to be tracked and accounted in iterative calculation. Therefore, 
the unlimited random structure shown in Fig{T]is not the actual lattice we are 
going to study, and two important principles have to be settled here: 1 ) the 
branch number on the vertices of one unit must be the same excluding the 
base vertex; 2 ) for any arbitrary square the three sub-trees contributing onto 
it towards to the original point should be identical. And the branch number 
on the vertices of different levels are random with the structural ratio P. A 
sample structure is presented in the Figj2^. 

Although these two limitations confine a locally ordered configuration on 
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Figure 2: The calculation scheme of the solution x and the free energy, (a) 
A demonstration of the actual structure employed in this paper to calcu¬ 
late the solutions. To ensure a local homogeneity to facilitate the recursive 
calculation, branch numbers on the three sites of one square are set to be 
identical. The recursive calculation towards an imaginary original point O. 
(b) The calculation of free energy by a replica counterpart to obtain smaller 
hooked-up lattices to yield the local free energy. 
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the same levels and impair the randomness of the lattice, from a general view 
of the infinitely large structure, we still have two or three branched vertex 
randomly appears with a fixed probability. Therefore, we may say that this 
paper discussed a special case of the ideally random multi-branched lattice 
with identical sub-trees contributions. 

The calculation method basically follows the previous works done by Gu- 
jrati’s 13, lla, [l6j, here we will briefly derive the method with specific modi¬ 


fications for the random-branched lattice. By defining the partition function 
of a sub-tree with the base spin fixed as +1 or —1 to be the partial partition 
function (PPF), denoted as Zf(± 1), where i is the level of the base spin, at 
the original site S$ (level 0 ) we have the total partition function: 


z» = z®(+ 1 ) + ZS(- 1 ), (5) 

The exponent B is the number of sub-trees (branches) contributing onto Sq. 
While in the regular lattice the branch number is a fixed variable, in this work 
the branch number is a series of random variables on each level, therefore 
the branch number shall be written as {B{\ associated to the level i. 

As shown in Fig J2^, for each square unit on level i it always takes in 
3(B i+ i — 1) sub-trees contributions from upper level and contribute itself to 
the unit on the lower level, then by defining B[ — B, — l, for the PPF on the 
level n we have the recursive relations 


z„<+) = ^z5r(s„, 2 )z5r(s„, s )zfi:r(s„,4)«>(7). (6) 

7=1 

16 

Zn(~) = J^Z*lt\S n i)Z*-i\S n3 )Z*;i\S nA )w{ 1 ). (7) 

7=9 

With the sites labeling shown in Figj2^,, on the two sites S n) 2 and S n $ neigh¬ 
boring to S n> 1 we have 2 B' n+2 sub-trees contributing from level n + 1 , and 
similar to the S n> 4 diagonal to S n . The w( 7 ) is the local weight of the square 
confined by the four sites Skp, S n> 2 , S n ^, and S nt 4 . In principle this local 
weight should exclude the weights of base site S n> i, e~ l3HSa ’ 1 , to conduct a 
cavity field contribution, however this can be ignored since we set the external 
held H = 0 in this paper. 

With the PPFs defined above, we introduce a ratio x{S n ) (simplified as 
x n ) on the site S n 
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( 8 ) 


_ Z n (+) 

Xn ~ z n (+) + z n (-y 

as the weights ratio of PPFs with the fixed spin state of S n , x denotes the 
sub-tree contribution to the magnetization held of the site S n , i.e. the cavity 
contribution (25|. Therefore this ratio can be treated as the “solution 1 ' of 
the system to indicate the magnetization of a site. Obviously, the solution 
of x = 0 or 1 indicates a uniform spins orientation to either up or down, and 
conversely the solution x = 0.5 implies a half-half probability to have a up or 
down spin on one site, i.e. the amorphous state without external magnetic 
held. 

Since x n is a function of PPFs on level n, and as shown in Eqj6]and[Zla 
PPF is a function of the PPFs on higher level, then x n can be derived as a 
function of xs on higher levels in a recursive fashion: 

X n = f{x n+ 1 ). (9) 

With an initial seeding input of x, the solution can be obtained by apply¬ 
ing the EqE] many times until reaching a fixed-point solution. Principally, 
various input seeds should be tried to exhaust all possible fixed points, how¬ 
ever this work deals with the simplest ferromagnetic Ising model, and a single 
solution at a fixed temperature is expected. 

The exact form of Eqj9] is derived below: Starting from 

Z n (+) 

Z n (-h) + Z n (—) 

Zn(~) 

Z n (+) + Z n (—) 


f x n if S n +1 
{ y n if S n = -1 


we define a compact note 


z n (S n ) 


A 


BL 


In terms of 


Z n (+) + Z n ( ), 





we have 


A 


B' 


»(=*=) = 


E 


D/ 

A ri+l 
/1 n+1 ' 


B n+1 

'n+1 


(S 




and 

*n(±) = ^^i 1 (^n 1 2 )^i 1 (5„ l3 )^t 1 (5 , „ I 4 )«;(7)/Q(a;n+i), 

where the sum is over 7 = 1,2 ,3,..., 8 for S Uj 1 = +1, and over 7 = 
9,10,11,..., 16 for = —1, and where 

Q( Xn+1 ) = a^/jCt 1 ; 
it is related to the polynomials 
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<ww = E 2 5t(s„,2)*5r(s„, 3 )*5r(s,, 4 )u>(7), 

7=1 

16 

Q-(x n+ 1 ) = ^^ 1 (5 , ni 2)^+t 1 (5 , ni 3)^t 1 (5 , n)4 )«;(7), 

7=9 


according to 

Q{^n+ 1) Q+{.^n+ 1) T Q— (Z/n+l) - 

In terms of the above polynomials, we can express the recursive relation for 
the ratio x n in terms of x n+ \: 


_ Q+(*^n+l) 

= Q(x n+1 ) ■ 

The above process is the static recursive calculation of solutions. For the 
random multi-branched lattice, on each iteration the program will randomly 
assign the value of B to be 2 or 3 according to the structural ratio P 2 , and then 
execute the calculation for that local level. Therefore, the solutions exhibit 
a fluctuation around an average value, and the fixed-point solution does 
not exist. According to our experience, the calculation reaches a “stable” 
solutions oscillation x n G (x ± a) in no more than 2000 iterations. In the 
program we did the calculation 16,000 times and average the last 11,000 xs 
to reach a reliable x. In this way, although the exact calculation is executed 
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in each iteration, instead of exact results what we actually obtained is a 
numerical solution distribution. 


2.3. Calculation of Thermodynamics 

The thermodynamic calculations again follows the same principle detailed 


in Ref. [13], 24], while a slight modification is necessary to fit the random 
conformation case. Herein we firstly review the general method to calculate 
the free energy on a homogeneous R-branched Husimi lattice. 

Imagine we achieved the fixed-point solution in the region around the 
original point O, the B squares joint on O are indexed as the first level, then 
there are 3 B x (B — 1) squares on the level 2. Now we cut off these sub-trees 
and hook them up to form 3 x (B — 1) smaller but identical lattices, the 
partition function of these lattices are 


Z 2 = Z"(+) + Z. 


The free energy of the left out local squares is 

F'local C Z^ B ' 

We have 4 /B cites in a square and B squares in the local origin region. The 
free energy per site is: 


4 ‘ 

By substituting Z n (+) = A n x n and Z n (—) = A n y n , we have 

1 0 2B ' 

F = ~A Tl ° g ^[ X B + (l- X )B]3B')- ( 10 ) 

An important concern in the inhomogeneous lattice discussed in this pa¬ 
per is that, the original point is merely a conceptual point that indicating the 
calculation direction, while it cannot be realized as an either 2 or 3 branched 
point. Therefore, the free energy of a local region should be calculated along 
the recursive process and averaged, in the same fashion of solutions, to ap¬ 
proach the free energy of the whole system. 

By this means, an replica counterpart calculation was designed to obtain 
the free energy along the recursive process, as shown in FigJ3). After the 
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solution converged to the stable oscillation, once an iteration of solution 
was calculated on one site, the immediate calculation of free energy of that 
site was carried out by hooking an identical replica half-tree onto that site 
(the shadow part in Figj2jo) . Then as in a complete lattice with the replica 
counterpart, the central point can be taken as the ‘origin’ and the general 
method of free energy calculation (Eq JTUD can be applied by cutting-off and 
rejoining the sub-trees on the corresponding sites S n and S' n to yield the free 
energy of the local region. Although this replica part does not exist in the 
system, due to the perfect symmetry it provides the correct free energy per 
site of the local region, and the value calculated on the real half is valid. 

We have also tested that, the branch number on the symmetric joint point 
does not affect the free energy calculation, i.e. hooking one or two imaginary 
counterparts onto the site to make R = 2 or 3 present the exact same free 
energy density, i.e. the thermodynamics only depends on the solution x once 
it was determined. 

3. Results and Discussion 

3.1. Spontaneous Magnetization 

The average solution x with P 2 = 0.5 is presented in FigJSk- The solution 
of deterministic cases of B = 2 and 3 with P 2 = 1 and 0 are also shown for 
comparison. As expected, the solution of P 2 = 0.5 with its deviation lies in¬ 
side the area lined out by the two deterministic solutions. Albeit the rigorous 
fluctuation in the intimidate temperature range, a distinct phase transitions 
can be figured out with the orientation preference of spins differentiated from 
0.5, i.e. the spontaneous magnetization of simple Ising model. 

Figure |3}o shows the detailed solution of P 2 = 0.5 around the critical tem¬ 
perature T c . Although the error bar is considerable large merely below T c , 
the transition is still distinctive. Besides the averaged solution over 11,000 
results, we have also investigated several single runs and found that the posi¬ 
tion of Tq is consistent, i.e. the phase transition of Ising spins is independent 
of the structural stochasticity. This phenomenon evidenced a well-defined 
thermodynamics of the randomly multi-branched lattice, and we can be con¬ 
fidential to apply the conventional analysis onto this inhomogeneous lattice. 

3.2. Thermal Behavior and Fluctuation 

Ten sets of free energies with P 2 varying from 2 to 3 by the increment 
of 0.1 are shown in FigJUu Similar to the solutions, free energies of cross- 
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0 1 2 3 4 5 


(a) 



3.35 3.40 3.45 3.50 3.55 3.60 3.65 3.70 3.75 3.80 


(b) 

Figure 3: (a) The solution x with P 2 = 1, 0.5 and 0. (b) The enlarged 
P 2 = 0.5 solution curve around the phase transition range. 
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dimensional structures lies inside the area outlined by the 2 and 3-branched 
free energy curves. The uneven curves of the free energies of models with 
2 < P < 3 imply the thermal fluctuations, which agree with the significant 
error bar of the solutions presented in Figj3l 

To more clearly investigate the thermal fluctuation, the entropies derived 
from S = —dF/dT of five structures are presented in FigHb. Rigorous fluc¬ 
tuations can be observed on the entropy behaviors of stochastic models, with 
a superheating extrapolation unphysically surpass the entropy of the amor¬ 
phous state. These singularities should not be understood as the simple effect 
of randomness, it actually reveals the inhomogeneity in the cross-dimensional 
range: For the near surface region of a bulk materials, the layer beneath the 
surface is at the halfway up the energy landscape, and its state is fragile and 
easily to be frustrated by the energy flow from surface to bulk or vice versa. 
At the critical range, the system is in dramatic dynamics, as the surface 
begins the order-disorder phase transition while the bulk remains ordered, 
therefore the energy fluctuations is expected to induce locally superheating 
in between. 

3.3. Branch Ratio vs. Transition Temperature 

One of the most important aim to design and investigate the inhomoge¬ 
neous Husirni lattice, is its application in the cross-dimensional range. Fig¬ 
ure [5^ presents the ten sets of solutions with various analog B, and the 3D 
mapping of these solutions is presented in FigJSJo. This 3D mapping clearly 
indicates the gradient in the cross dimensional situation, i.e. the depth from 
surface to the bulk, or the expansion of a thin film to thick bulk. Herein we 
can observe from the 3D mapping that it is possible to have different phases 
on the near-surface region at the same temperature. The edge of the 0.5 
plateau is the transition line, which is fitted as 

T c = -0.23436 + 1.04126e°' 53187B (11) 

with R 2 = 0.99991 (FigJSfc). 

Before analyzing the exponential behavior of T c vs. B presented in EqJTTl 
the variation of ground state energy with the analog B should be figured 
out as a reference. When the systems is at T = 0, all the spins point 
to one direction as a perfect crystal and we have the ground state energy 
E = F with entropy S — 0. For the homogeneous lattices of B = 2 and 
3, the ground energy state is simply the number of coupling interactions 
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(a) 



(b) 

Figure 4: (Color online) (a) The free energy of models with analog B varies 
from 2 to 3. (b) The entropy of models with B — 2, 2.2, 2.5, 2.7 and 3. 
The two deterministic models presents a smooth boundary line, while the 
stochastic models exhibit rigorous thermal fluctuations in between. Clear 
singularities can be observed on these curves around the critical tempera¬ 
ture, implies the inhomogeneous phases in the cross-dimensional range, where 
superheating are possible. 
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Figure 5: (Color online) The transition temperature Tc variation with branch 
ratio, the ratio is presented by the analog branch number B by Eq.l: (a) ten 
sets of solutions of B from 2 to 3 with 0.1 increment; (b) the 3D mapping 
of the ten sets of solutions in the T-X-B space; (c) the transition line as the 
edge of the 0.5 plateau extracted from the 3D mapping. 

15 













Table 1: Variation of T c and the ground state energy density with analog B. 


analog B 

T 

O 

II 

2 

2.78 

-2 

2.1 

2.95 

-2.103 

2.2 

3.12 

-2.203 

2.3 

3.31 

-2.301 

2.4 

3.49 

-2.410 

2.5 

3.70 

-2.506 

2.6 

3.92 

-2.606 

2.7 

4.15 

-2.703 

2.8 

4.37 

-2.802 

2.9 

4.64 

-2.902 

3 

4.90 

-3 


Et=o,b= 2 = 4 x y = —2, or Et=o,b=2 = 6 x y = -3 (one interaction 
is shared by two spins). For the cross-dimensional lattices in between, the 
ground state energies Et =o with various analog B exhibit linear behavior 
as presented in Tabled! this phenomenon evidenced that the system has a 
well-defined average property according to the probability ratio. 

Therefore, the exponential correlation of T c vs. analog B suggests that 
the stochasticity might be the single reason to affect the transition behav¬ 
ior near the critical temperature, and we may hypothesize that structural 
randomness reduces the stability of the system and leads to a lower transi¬ 
tion temperature. Nevertheless, in the quantitative aspect, the exponential 
fitting and the meaning of the parameters in Eq{lT] is still unclear. Further 
investigations on this phenomenon is in progress. 

4. Conclusion 

The ferromagnetic Ising model was solved on a randomly multi-branched 
Husirni square lattice. The lattice was constructed with two or three squares 
connecting on one with the probability P 2 or P 3 . While the system is de¬ 
terministic on each end, the structures with variable P in between describes 
a cross-dimensional situation from 2D to 3D, such as the thin film or the 
near-surface case. The solution x representing the magnetization of spins 
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was obtained by the cavity partial partition function method, and the ther¬ 
modynamics such as free energy and entropy were derived from the solutions. 

Typical spontaneous magnetization were observed on the random inho¬ 
mogeneous lattices, the transition occurs on an characteristic Tq regardless 
of the randomness and thermal fluctuation. The Tq , however, behaves an 
exponential correlation to the structural ratio or analog branch number. Con¬ 
sidering the linear relationship of the ground state energy E ?=o vs. B, the 
cause of these lowered X^s is hypothesized to be the stochasticity that re¬ 
duced the stability of the system, while the quantitative meaning of the 
fitting equation remains unclear. Rigorous fluctuations were observed on the 
entropy behaviors with singularities of locally supercooling or superheating 
extrapolation, implying the severe inhomogeneity in the cross-dimensional 
range around the critical temperature. 
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